#ifndef am_InterVar_H
#define am_InterVar_H

#include <string>
#include <iostream>
#include <tbb\concurrent_vector.h>
#include <boost\shared_ptr.hpp>

#include <borealis\Raster.hpp>
#include <borealis\AsciiRaster.hpp>
#include <borealis\Interpolation.hpp>
#include "..\AuroraTypes.hpp"

namespace Aurora
{
	/*
	Class: am_InterpVar
	Interpolates a variable over the domain

	In Variables:
	Raster<float>	DEM variable
	string			Varaible to interpolate
	string			Method
						- "IDW"
						- "TPST"
	Out Variables:
	string			Name of gridded variable


	*/
	class am_InterpVar : public AuroraTypes::ModuleBase
	{
	public:

		//in variables:
		//DEM
		//var
		//interp method (IDW is only supported)
		int DEM;
		int VAR;
		int METHOD;

		//out variable
		//Interpolated  raster
		int OUTVAR;

		am_InterpVar(std::string n)
		{
			DEM			=0;
			VAR			=1;
			METHOD		=2;

			OUTVAR		=0;

			name = n;
		}



		void run(AuroraTypes::StationList* stations,AuroraTypes::VariableList* variables, AuroraTypes::InVariables& inVars, AuroraTypes::OutVariables& outVars)
		{

			boost::shared_ptr<Borealis::Raster<float> > dem;
			boost::shared_ptr<Borealis::AsciiRaster> interp;

			//bail if we don't have any met stations
			if(stations->size() == 0)
				throw std::runtime_error("In " + std::string(__FILE__) +" @ line " + boost::lexical_cast<std::string>(__LINE__) +"No met data");

			try
			{
				dem = GetRaster(variables,inVars[DEM]);

			}
			catch (std::runtime_error e)
			{
				throw std::runtime_error("In " + std::string(__FILE__) +" @ line " + boost::lexical_cast<std::string>(__LINE__) +" - failed to get DEM parameter");
			}


			interp.reset(new Borealis::AsciiRaster(dem->GetRow(),dem->GetCol()));
			interp->SetXllCorner(dem->GetXllCorner());
			interp->SetYllCorner(dem->GetYllCorner());
			interp->SetCellSize(dem->GetCellSize());


			//for each station lower the temperature to a common level z
			//build a new point and save it for the interpolation scheme
			tbb::concurrent_vector<Borealis::RasterPoint<float> > points;

			//correct the reference temperatures
			for(unsigned int i = 0; i < stations->size();i++)
			{
				double variable = stations->at(i)->Now().Get<float>(inVars[VAR]);
				unsigned int col = stations->at(i)->GetX();//col
				unsigned int row = stations->at(i)->GetY();//row

				Borealis::RasterPoint<float> point(row,col,variable);
				points.push_back( point );
			}

			std::string interpMethod = inVars[METHOD];
			if(interpMethod == "IDW")
			{
				Borealis::Interpolation::InverseDistanceWeighting(points, interp);
			}
			else if(interpMethod ==  "TPST")
			{
				Borealis::Interpolation::ThinPlateSplineT(points,interp);
			}
			else
			{
				throw std::runtime_error("Unknown interpolation method");
			}

			AuroraTypes::VariableList::accessor a2;
			variables->insert(a2,outVars[OUTVAR]);		
			a2->second = interp;

		}

	};
}

#endif